;; Non-robust Graham Scan algorithm for computing the convex hull of a set of 2D points.
;; Copyright (C) 2021  Rocks Mazama
;;
;; This program is free software: you can redistribute it and/or modify
;; it under the terms of the GNU General Public License as published by
;; the Free Software Foundation, either version 3 of the License, or
;; (at your option) any later version.
;;
;; This program is distributed in the hope that it will be useful,
;; but WITHOUT ANY WARRANTY; without even the implied warranty of
;; MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
;; GNU General Public License for more details.
;;
;; You should have received a copy of the GNU General Public License
;; along with this program.  If not, see <https://www.gnu.org/licenses/>.

(ns mazama.compgeom.convex-hull.graham-scan
  (:require [mazama.compgeom.core :refer [distance-squared point-orient2d]]))

;;;; NOTE
;;;;
;;;; 1. This software may fail if floating point values are used.
;;;;
;;;; 2. This software assumes that:
;;;;    the x-axis increases to the right;
;;;;    the y-axis increases DOWNWARD.

(defn- bottom-left-point-fn [[x1 y1 :as p1] [x2 y2 :as p2]]
  (cond
    (> y1 y2) p1
    (> y2 y1) p2
    (< x1 x2) p1
    :else p2))

(defn- bottom-left-point [S]
  (reduce bottom-left-point-fn S))

(defn- compare-fn
  "Point `b` is less than point `c` if `b` is further clockwise from `a` than `c` is."
  [a b c]
  (cond
    (= a b) -1
    (= a c) 1
    (= b c) 0
    :else (let [status (point-orient2d a b c)]
            (if (== status 0)
              (let [dist-ab (distance-squared a b)
                    dist-ac (distance-squared a c)]
                (if (< dist-ab dist-ac) -1 1))
              (- status)))))

(defn- reduce-fn
  "Add point `c` to the `hull`, discarding previous points on the hull if necessary."
  [hull c]
  (if (< (count hull) 2)
    (conj hull c)
    (let [hull- (pop hull)
          b (peek hull)
          a (peek hull-)]
      (case (point-orient2d a b c)
        -1 (recur hull- c)
        0 (recur hull- c)
        1 (conj hull c)))))

(defn convex-hull
  "Compute the convex hull of the 2D points in `S`.
  https://en.wikipedia.org/wiki/Graham_scan"
  [S]
  (if (< (count S) 3)
    S
    (let [initial-point-on-hull (bottom-left-point S)
          P (sort (partial compare-fn initial-point-on-hull) S)]
      (reduce reduce-fn [] P))))
